#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Sun Jun 30 22:56:11 2024

@author: liqingsimac
"""
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
fig=plt.figure()
ax=plt.axes(projection='3d')
x=np.linspace(0,1,200)
theta=np.linspace(-2*np.pi,15*np.pi,300)
[x,theta]=np.meshgrid(x,theta)
phi=(np.pi/2)*np.exp(-theta/8/np.pi)
X=1-.5*pow(1.25*pow(1-np.mod(3.6*theta,2*np.pi)/np.pi,2)-1/4,2)
y=1.95653*pow(x,2)*pow(1.27689*x-1,2)*np.sin(phi)
r=X*(x*np.sin(phi)+y*np.cos(phi))
h=X*(x*np.cos(phi)-y*np.sin(phi))
#c=cm.get_cmap('magma')
c=cm.get_cmap('winter')
surf=ax.plot_surface(r*np.cos(theta),r*np.sin(theta),h,rstride=1,cstride=1,cmap=c,linewidth=0,antialiased=True)
ax.set_position((0,0,1,1))

